Purpose

Create maps to visualize median household income across New York City using ACS data.

Set-up

# Enable caching of data
options(tigris_use_cache = TRUE)

Packages

# Loading packages
library(rmarkdown) # markdown package to help knit to git
library(dplyr)   # for data manipulation
library(purrr)  # map list & other functions
library(ggplot2) # for graphing
library(sf)     # simple features - encoding spatial vector data
library(svglite) # for exporting SVG files
library(stringr) # text/string functions
library(shiny) # interactive dashboard
library(tidycensus) # tidyverse census functions
library(tigris) # read census shapefiles

library(leaflet) # interactive maps
library(ggiraph) # interactive maps via ggplot2

Census API key

In order to access the Census API, you need to enter a Census API key, which you can request here: https://api.census.gov/data/key_signup.html

# Use function to install your API key for future sessions. This will add your CENSUS API key to your .Renviron file so it can be called securely without being stored in your code. 
census_api_key("[INSERT YOUR KEY]", install = TRUE)

# Alternatively, access the API by entering in your API key in quotes below just once each session
Sys.setenv(CENSUS_API_KEY = "[INSERT YOUR KEY]")
# After you have installed your key, it can be called any time by typing Sys.getenv("CENSUS_API_KEY") and can be used in package functions by simply typing CENSUS_API_KEY
Sys.getenv("CENSUS_API_KEY")

Load ACS Data

Will use the 5-year estimate data from 2015 for median household income to start.

Output variable

# Load variables from ACS5 2015 census and isolate the variable for median household income
vacs15 <- load_variables(2015, "acs5", cache = T)

# Confirm this is the correct outcome variable needed
vacs15 %>% 
  filter(str_detect(name, "B19113_001"))
## # A tibble: 1 × 4
##   name       label                                             concept geography
##   <chr>      <chr>                                             <chr>   <chr>    
## 1 B19113_001 Estimate!!Median family income in the past 12 mo… MEDIAN… tract
# Good - label confirms this is the correct variable    

Income from ACS5 by Census Tract Table

Load income data by tract and subset to the 5 NYC boroughs.

incmedhh_2015 <- get_acs(geography = "tract",
                             variables = c(incmedhh = "B19113_001"), 
                             state = "NY",
                             county = c("New York", 
                                        "Kings", 
                                        "Queens", 
                                        "Bronx", 
                                        "Richmond County"),
                             year = 2015,
                             geometry = T)

# Preview df
head(incmedhh_2015)
## Simple feature collection with 6 features and 5 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -73.92922 ymin: 40.80354 xmax: -73.87251 ymax: 40.8299
## Geodetic CRS:  NAD83
##         GEOID                                       NAME variable estimate  moe
## 1 36005002300    Census Tract 23, Bronx County, New York incmedhh    17671 6565
## 2 36005002701 Census Tract 27.01, Bronx County, New York incmedhh    24784 6865
## 3 36005004100    Census Tract 41, Bronx County, New York incmedhh    18671 6541
## 4 36005004800    Census Tract 48, Bronx County, New York incmedhh    31538 6814
## 5 36005006700    Census Tract 67, Bronx County, New York incmedhh    21747 5670
## 6 36005008900    Census Tract 89, Bronx County, New York incmedhh    30385 5751
##                         geometry
## 1 MULTIPOLYGON (((-73.92922 4...
## 2 MULTIPOLYGON (((-73.91796 4...
## 3 MULTIPOLYGON (((-73.92468 4...
## 4 MULTIPOLYGON (((-73.8764 40...
## 5 MULTIPOLYGON (((-73.92178 4...
## 6 MULTIPOLYGON (((-73.89704 4...

Static map with ggplot2

Get shapefile for NYC water layer. Will use as separate layer for map formatting.

nyc_water <- area_water("NY", c("New York", 
                                "Kings", 
                                "Queens", 
                                "Bronx", 
                                "Richmond County"))

Plot data using ggplot2

ggplot(incmedhh_2015) +
  # Income data
  geom_sf(alpha = 1, 
          color = NA, 
          aes(fill = estimate)) +
  # Theme removes cartographic background grid & axis labels
  theme_void() +
  # Legend
  scale_fill_viridis_c(labels = scales::label_dollar(),
                       name = "Median Household Income") + 
  # Remove water layer by setting color to white
  geom_sf(data = nyc_water, 
          fill = "white",
          color = NA, 
          aes(geometry = geometry)) +
  # Set background to white to match water layer
  theme(panel.background = element_rect(fill = "white", 
                                        color = NA)) 

Interactive maps

Create interactive maps using leaflet and ggiraph packages

leaflet

# Set color palette argument
pal <- colorNumeric(
  palette = "magma",
  domain = incmedhh_2015$estimate
)

# Plot
leaflet() %>% 
  #providers argument sets map background. Chose one of the free versions available through leaflet
  addProviderTiles(providers$CartoDB.Positron) %>%   
  addPolygons(data = incmedhh_2015,
              color = ~pal(estimate),
              weight = 0.5,
              smoothFactor = 0.2,
              fillOpacity = 0.5,
              label = ~estimate)

ggiraph

ggiraph is a package that makes ggplot2 figures interactve. Use the basic code set up from the original static ggplot2 map for this section and edit as needed.

# Create new df that has a text column combining each tract name and the income amount (i.e. "[TRACT NAME]: [$ AMOUNT]"). This will appear when cursor hovers over each tract
incmedhh_2015_ggir <- incmedhh_2015 %>% 
  mutate(tooltip = paste(NAME, estimate, sep = ": "))

# Plot with ggiraph
ggir <- ggplot(incmedhh_2015_ggir, aes(fill = estimate)) +
  
  # Follow same set up as ggplot2 static map, but with interactive version of ggplot2 functions
  geom_sf_interactive(aes(tooltip = tooltip, data_id = NAME), 
                      size = 0.1) + 
  scale_fill_viridis_c_interactive(option = "plasma",
                                   labels = scales::label_dollar()) +
  # Color water layer white
  geom_sf(data = nyc_water, 
          fill = "white",
          color = NA, 
          aes(geometry = geometry)) +
  labs(title = "Median Household Income by Census Tract",
       caption = "Data source: 2015 5-year ACS, US Census Bureau",
       fill = "ACS estimate") +
  # Theme removes cartographic background grid & axis labels
  theme_void() +
  # Set background white to match water layer
  theme(panel.background = element_rect(fill = "white", 
                                        color = NA)) 

# Output
girafe(ggobj = ggir) %>% 
  # Fill tract with select color when cursor hovers over it
  girafe_options(opts_hover(css = "fill:cyan;"),
                 opts_zoom(max = 10))

Multi-year interactive map

Create interactive map that uses data from multiple years from the ACS-5 survey

# Save years of interest to list
years <- lst(2010, 2015, 2020)

# Run function to load census data across years in list
incmedhh_all <- map_dfr(
  years,
  ~ get_acs(
    geography = "tract",
    variables = "B19113_001",
    state = "NY",
    county = c("New York", 
               "Kings", 
               "Queens", 
               "Bronx", 
               "Richmond County"),
    year = .x,
    survey = "acs5",
    geometry = T
    ),
  .id = "year"
  ) %>% 
  mutate(year = as.numeric(year))

Shiny

Create interactive dashboard with Shiny to look at multi-year data using ggiraph package

# Create text column for hover over feature
incmedhh_all_ggir <- incmedhh_all %>% 
  mutate(tooltip = paste(NAME, estimate, sep = ": "))
# Shiny setup

ui <- fluidPage(
  # Create sidebar with dropdown menu to select year
  sidebarLayout(
    sidebarPanel = sidebarPanel(
      selectInput(
        'selected_year',
        'Select a year to map',
        choices = unique(incmedhh_all_ggir$year)
      )
    ),
    mainPanel = mainPanel(
      girafeOutput('girafe_output', height = 600)
    )
  )
)
      
server <- function(input, output, session) {
  
  # Reactive function that filters by year based on input
  dat_year <- reactive({incmedhh_all_ggir |> 
      filter(year == input$selected_year)})
  
  output$girafe_output <- renderGirafe({
    
    # ggiraph/ggplot2 code (same setup as the previous sections)
    ggir_all <- ggplot(data = dat_year(), 
                       aes(fill = estimate)) +
      geom_sf_interactive(aes(tooltip = tooltip, data_id = NAME), 
                         size = 0.1) +
      scale_fill_viridis_c_interactive(option = "plasma",
                                       labels = scales::label_dollar()) +
      # Color water layer white
      geom_sf(data = nyc_water, 
              fill = "white",
              color = NA, 
              aes(geometry = geometry)) +
      labs(title = "Median Household Income by Census Tract",
           caption = "Data source: 5-year ACS, US Census Bureau",
           fill = "ACS estimate") +
      # Theme removes cartographic background grid & axis labels
      theme_void() + 
        # Set background white to match water layer
      theme(panel.background = element_rect(fill = "white", 
                                            color = NA)) 

    girafe(ggobj = ggir_all,
           options = list(opts_hover(css = "fill:cyan"),
                          opts_zoom(max = 10)))
      
    
  })
  
  # Shiny session runs until you end it by closing pop up window
  session$onSessionEnded(function() {
    stopApp()
  })
  
} 

shinyApp(ui = ui, server = server)
Shiny applications not supported in static R Markdown documents